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! Abstract 

In the present article, we present a new gravitational galactic model, describing motion in elliptical as well as in disk 
galaxies, by suitably choosing the dynamical parameters. Moreover, a new dynamical parameter, the S{g) spectrum, is 
introduced and used, in order to detect islandic motion of resonant orbits and the evolution of the sticky regions. We 
investigate the regular or chaotic character of motion, with emphasis in the different dynamical models and make an 
extensive study of the sticky regions of the system. We use the classical method of the Poincare (r — p r ) phase plane 
and the new dynamical parameter of the S(g) spectrum. The L.C.E is used, in order to make an estimation of the 
'Hepree of chaos in our galactic model. In both cases, the numerical calculations, suggest that our new model, displays 
a wide variety of families of regular orbits, compared to other galactic models. In addition to the regular motion, this 
rHi. new model displays also chaotic regions. Furthermore, the extent of the chaotic regions increases, as the value of the 
flatness parameter b of the model increases. Moreover, our simulations indicate, that the degree of chaos in elliptical 
i-C galaxies, is much smaller than that in dense disk galaxies. In both cases numerical calculations show, that the degree 
chaos increases linearly, as the flatness parameter b increases. In addition, a linear relationship between the critical 
Q value of angular momentum and the b parameter if found, in both cases (elliptical and disk galaxies). Some theoretical 
arguments to support the numerical outcomes are presented. Comparison with earlier work is also made. 

& "Keywords: Galaxies: kinematics and dynamics; new dynamical models 
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1. Introduction 

It is very interesting, that scientists interested in Dy- 
namical Astronomy, are often model constructors. Dy- 
namical models, are very useful tools for the study of the 
behavior of orbits in a stellar dynamical system, such as a 
galaxy. A dynamical model is usually a mathematical ex- 
pression, giving the potential as a function of the distance 
from center of the galaxy. Although, sometimes instead of 
the potential, the density or the gravitational force as a 
function of the distance, is given. 

Important and also interesting models, have been built 
up during the last three decades, when the idea that the 
elliptical galaxies were considered as rotationally flattened 
systems, has been abandoned. Two main categories of 
models, were used for the study of the galactic motion. 
(1) The local models, built up of harmonic oscillators (see 
Saito and Ichimura 1974, Caranicolas 1984, 1993, 1994, 
Deprit 1991, Elipe et al. 1995, Elipe 1999, Karanis and 
Caranicolas 2002) and (2) The global mass models, de- 
scribing global motion (see Clutton-Brock et al. 1976, 
Martinet and Pfenniger 1986, Carlberg and Innanen 1987, 



'Corresponding author: 
E-mail address: evzotos@astro.auth.gr (Euaggelos E. Zotos) 



Patsis and Zachilas 1990, Flynn et al. 1996, Fellhauer et 
al. 2006). A detailed description for both the above two 
kinds of galactic dynamical models can been found, by the 
reader in Binney and Tremaine (2008). 

In several papers, we see that low angular momen- 
tum stars in galaxies with massive nuclei, are scattered to 
much higher scale heights, displaying chaotic motion (see 
Caranicolas & Innanen 1991, Caranicolas 1997, Karanis 
& Caranicolas 2001, Caranicolas & Papadopoulos 2003, 
Caranicolas & Zotos 2010). All the above research, not 
only describes the chaotic character of orbits of low angular 
momentum stars, but also gives the main reason respon- 
sible for this behavior. The reason is the presence of the 
strong force, which comes from the massive nucleus. The 
dynamical models, used in all previous papers were global 
models. In Caranicolas & Innanen (1991) and Caranicolas 
(1997), they have used the mass model introduced by Carl- 
berg & Innanen (1987), in Karanis & Caranicolas (2001), 
they have used a logarithmic potential, while in Caranico- 
las & Papadopoulos (2003) a time - dependent model was 
used. 

The most characteristic result, derived in all cases, was 
the linear dependence of the critical value of angular mo- 
mentum and the mass of the nucleus. From our point of 
view, the background of this series of paper, in combina- 
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Figure 2: (a-b): (a-left): Density contours for the elliptical galaxy and (b-right): Density contours for the disk galaxy. The values of all 
parameters arc given in text. 
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Figure 3: (a-b): The mass density in the galactic plane z = 0, (a- left): for the elliptical galaxy and (b-right) for the disk galaxy. 



tion of other orbit calculations, was that managed to find namical system, such as the angular momentum and the 
and present, both numerically and semi - theoretically, re- degree of chaos in elliptical and disk galaxies, 
lationships connecting the physical parameters of the dy- The paper is organized as follows: In Section 2, we 
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present our new dynamical model, which is a combination 
of the logarithmic model and the Miyamoto - Nagai model 
(1975). In Section 3, we study the character of orbits 
in two different cases (i) Our model describes an ellipti- 
cal galaxy and (ii) our model describes a flat disk galaxy. 
Moreover, a comparison between those two cases is inves- 
tigated. In Section 4, we introduce and use the new S(g) 
spectrum. In the next Section 5, we present some semi- 
theoretical arguments, in order to support our numerically 
obtained outcomes. Finally, in the last Section a discus- 
sion and the conclusions of this research are presented. 

2. Presentation of the dynamical model 

In the present paper, we shall study the properties of 
motion in a galaxy described, by the potential 



V g {r,z) 



■In 



r 2 + 



(a+y/h? + bz 2 y 



(1) 



were r, z are the cylindrical coordinates, c is the scale 
length of the galaxy, while a and ft are parameters regard- 
ing the structure of the dynamical system. The parameter 
wo, is used for the consistency of the galactic units. The 
parameter b is very important, because it represents the 
flatness of the galaxy. The shape of galaxy is spherical or 
oblate when 1 < b < 2 or prolate when 0.1 < b < 1. In 
the present paper, we shall study in detail, the case of two 
oblate galactic systems (disk and elliptical galaxy) with 
the same value of flatness parameter (b = 1.3). Moreover, 
we shall compare the relationships between other parame- 
ters of the dynamical systems and the flatness parameter, 
when 1 < b < 2. There are three main line of arguments, 
justify our choice: (i). It is a simple dynamical model, 
describing motion in a galaxy with an axis and a plane 
of symmetry, (ii). This potential can describe motion, in 
multiple models, from an elliptical galaxy (a -C b), to a flat 
disk system, (hi). Combining potential (1) with the po- 
tential of a spherically symmetric nucleus, we can obtain 
a new model, displaying a wide variety of regular orbits, 
together with chaotic regions. 

Adding to potential (1) a spherically symmetric nu- 
cleus, we obtain the total potential 



Vtot (r, z) = V g (r, z) + V n (r, z) , 



where 



V n (r, z) = 



yV 2 + z 2 + c 2 



(2) 



(3) 



while, M n is the mass and c„ represents the scale length 
of the nucleus. The Plummer sphere we use in order to 
increase the central mass, has been applied many times in 
the past, in order to study the effect of the introduction 
of a central mass component in a galaxy (see Hasan and 
Norman 1990, Hasan et al. 1993). In this work, we use a 
system of galactic units, where the unit of mass is 2.325 x 
10 7 M©, the unit of length is lkpc and the unit of time is 



0.997748 x 10 8 yr. The velocity unit and the energy unit 
(per unit mass) are 10km/ s and 100(fcm/s) 2 respectively, 
while G is equal to unity. In these units, we use the values: 
c = 2.5, M n = 400, c„ = 0.25, while vo,a,h 7 b are treated 
as parameters. 

Figure la, shows the rotation curve of the galaxy de- 
scribed by potential (2), when: v n = 15, a = 0.1, b = 
1.3, h = 0.12. In this case, we have an elliptical galaxy, 
with a small circular velocity, except in the central region 
and a low mass density distribution. At r = 8.5, z = 0, 
we find a circular velocity 9 = 160km /s and a mass 
density p = O.O14M /pc 3 . On the contrary, when vq = 
25, a = 1.2,6 = 1.3, ft = 0.2, equation (2) describes a flat 
disk system. The corresponding rotation curve is shown 
in Fig. lb. Here things are very different, because at 
r = 8.5, z = 0, we obtain a circular velocity O = 246km/ s 
and a mass density p — O.134M0/pc 3 . It is very useful, 
to compute the galactic density p{r,z), derived from the 
total potential (2) using the Poisson's equation 



P {r, z) 



1 



4ttG 



W tot (r,z). 



(4) 



Figures 2a-b, shows the contours p (r, z) — const, for 
the elliptical and disk galaxy respectively. For the ellipti- 
cal galaxy in Fig. 2a, the contours are: (0.057, 0.021, 0.01, 
0.0058, 0.0037, 0.0025), while for the disk model in Fig.2b 
the contours are: (0.34, 0.14, 0.076, 0.044, 0.028, 0.02). In 
Figures 3a-b, we can observe the plots of the mass density 
in the galactic plane p (r, z — 0), as a function of the radius 
r, for the elliptical and disk galaxy respectively. Due to the 
fact that, the potential is logarithmic, for large values of r 
and z, the mass density varies like 1/r 2 and 1/z 2 respec- 
tively. This means that the total mass M(R), enclosed 
in a sphere of radius R increases linearly with distance. 
This explains why circular velocity profiles, shown in Fig. 
la-b are flat. Here we must point out, that our gravita- 
tional potential is truncated at R max — 20kpc, otherwise 
the total mass of the galaxies modeled by this potential 
are infinite, which is obviously not physical. 

3. Structure of the phase plane - Orbit calculations 

As the total potential V to t {r, z) is axially symmetric 
and the L z component of the angular momentum is con- 
served, the dynamical structure of the system, can be stud- 
ied using the effective potential 



Veff (r, z) = ^§ + V^t (r, z) , 



(5) 



in order to study the motion in the (r, z) plane. The equa- 
tions of motion are 



r=p r , z=p z 



Pr 



dv 



eff 



n 7 Pz — ! 

or oz 
where the dot indicates derivative with respect to time. 



d V 



eff 



(6) 
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The corresponding Hamiltonian can be written as 

H=±(p 2 r +p 2 z )+V eff (r,z) = E, (7) 

where p r and p z , are the momenta per unit mass, conjugate 
to r and z respectively, while E is the numerical value of 
the Hamiltonian, which is conserved. Equation (7) is an 
integral of motion, which indicates that the total energy of 
the test particle is conserved. The Hamiltonian (7), also 
describes the motion in the (r, z) meridian plane, rotating 
at the angular velocity 

* = « = (8) 

In this Section, we shall use the classical method of 
the Poincare (r,p r ), z = 0,p z > phase plane, in order to 
determine the regular or chaotic nature of motion. If we 
set z = p z = in equation (7), we obtain the limiting curve 
in the (r — p r ) phase plane, which is the curve containing 
all the invariant curves, for a given value of the energy 
integral E. The limiting curve of the dynamical system 
corresponds to 

\p 2 r + V ef f{r)=E. (9) 

Orbit calculations, are mainly based on the numerical 
integration of the equations of motion (6), which was made 
using a Bulirsh - Stoer numerical integration Fortran code, 
with double precision subroutines. The accuracy of the 
calculations, was always checked by the constancy of the 
energy integral, which was conserved up to the twelfth 
significant figure. 

Two cases will be studied in this investigation: (a) The 
model (1) represents an elliptical galaxy and (b) The 
model (1) represents a flat disk galaxy. 

3.1. Elliptical galaxy model 

Figure 4 shows the structure of the Poincare phase 
plane, when v = 15, a = 0.1,6 = 1.3, ft = 0.12, for the 
Hamiltonian (7). It is evident that, for the above values of 
the dynamical parameters, we study the orbital structure 
of an elliptical galaxy. The value of the energy is E = 486 
and the value of the angular momentum is L z = 10. 

One can observe large areas of regular motion and a 
unified chaotic sea. There are three distinct regions of 
regular orbits, (i) Orbits forming invariant curves sur- 
rounding the central invariant point, (ii) orbits forming 
the triple set of islands and (hi) orbits forming the dou- 
ble set of islands symmetric with respect to the r - axis. 
The corresponding orbits are shown in Figure 5a-d. Fig 5a 
shows an orbit forming an island near the central invariant 
point. The initial conditions are: r = 2.8, p r o = 0, z = 0, 
while for all orbits the value of p z o, is always found using 
the energy integral (7). This is a typical tube orbit. The 
orbit shown in Fig. 5b, belongs to group (ii) and its initial 
conditions arc: r = 9.8, p r o = 0,z = 0. This orbit is 
a characteristic example of the 4/3 resonance. The orbit 
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Figure 4: The (r — p r ) phase plane for the elliptical galaxy, when 
E = 486. Sec text for details. 




2 4 6 8 10 



Figure 6: The (r — p r ) phase plane for the disk galaxy, when E = 
1425. See text for details. 

shown in Fig. 5c belongs to group (hi). The initial condi- 
tions of this orbit are: ro = 1.4, p r o = 28.5, zq = 0. This is 
a banana orbit and carries stars out of the galactic plane. 
Finally, the orbit shown in Fig. 5d is a chaotic orbit with 
initial conditions: ro = 2.3, p r o = 24, zq = 0. All orbits 
were calculated for a time period of 100 time units. 

3.2. Disk galaxy model 

Figure 6 shows the structure of the corresponding phase 
plane when vq = 25, a = 1.2,6 = 1.3, ft = 0.2. This time 
we have a disk galaxy. The value of the energy is E = 1425 
and the value of the angular momentum is L z = 10. 

Once more, one can see areas of regular motion and a 
large chaotic sea. The regular orbits belong to four main 
families, (i) Orbits forming invariant curves surrounding 
the central invariant point, (ii) orbits forming the four is- 
lands symmetric to the r - axis, (hi) orbits forming the 
triple set of islands and (iv) orbits forming the outer dou- 
ble set of islands. The corresponding orbits are show in 
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Figure 5: (a-d): Orbits in the elliptical galaxy model, when E = 486. The values of all parameters are as in Fig. 4. 
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Figure 7: (a-d): Orbits in the disk galaxy model, when E = 1425. The values of all parameters are as in Fig. 6. 



Figure 7a-d. Figure 7a shows an orbit forming an island 
near the central invariant point. The initial conditions are: 
ro = 1.7, Pro = 0, zo = 0. In Fig. 7b we can see a resonant 
orbit, with initial conditions: ro = 2.8, p r o = 0, zq = 0. 
The orbit in Fig. 7c belongs to group (ii). The initial 
conditions are: ro = 1.2, p r o = 11, zq = 0. The orbit 
shown in Fig. 7d is a chaotic orbit with initial conditions: 
r = 0.4, p r0 = 0,z = 0. All orbits were calculated for a 



time period of 100 time units. 

3. 3. Comparison of the two models 

It is important to note that in the case of the disk 
galaxy, the phase plane of the dynamical system displays 
a wider variety of regular orbits, than in the case of the 
elliptical galaxy. Moreover, the extent of the chaotic sea 
in the phase plane of the disk galaxy model (Fig. 6) is 
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Figure 8: The L.C.Es for a chaotic orbits in elliptical and disk galaxy 
models. 
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Figure 9: Maximal L.C.Es of a chaotic orbit in elliptical and disk 
galaxy vs. b. 
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Figure 10: Plot of the percentage of the area A% in the phase plane 
covered by chaotic orbits vs. b for the elliptical and disk galaxy 
models. 



significant larger than this of the elliptical galaxy (Fig. 4) . 
In order to have a better estimation of the degree of chaos, 
we have computed the maximal Lyapunov Characteristic 
Exponent (L.C.E) (see Lichtenberg and Lieberman, 1992 
for details), for a chaotic orbit in each case. The results 
are shown in Figure 8. The initial conditions of the chaotic 
orbit for the elliptical galaxy are: r = 2.5, p r o = 21, Zo = 
0, while for the disk galaxy the chaotic orbit has initial 
conditions: tq = 9.7, p r o = 0,zq = 0. Both orbits were 
calculated for a time period of 10 6 time units. The values 



Figure 12: Plot of the L zc vs. b, for both dynamical models. 

of all parameters are as in Fig. 4 (for the elliptical galaxy) 
and Fig. 6 (for the disk galaxy). It is evident, that the 
L.C.E of the disk system is twice larger than the L.C.E 
of the elliptical galaxy. Therefore, one may conclude, that 
the degree of chaos in disk galaxies is larger, compared to 
that observed in elliptical galaxies. 

Let us now come to see how the flatness parameter b of 
our model, is connected with other important parameters 
of the dynamical system. In all cases the value of the flat- 
ness parameter b, is between 1 < b < 2. Figure 9 shows a 
plot of the value of the maximal Lyapunov Characteristic 
Exponent (L.C.E) of a chaotic orbit in elliptical and disk 
galaxy vs b. One can see, that the values of the L.C.E 
increases linearly as b increases in both cases. Figure 10 
shows the relationship between the percentage of the area 
A% covered by chaotic orbits in the phase plane and the 
flatness parameter b. We observe, that in both cases (ellip- 
tical and disk galaxy) the percentage A% increases as the 
flatness parameter b increases. Since A% cannot exceed 
100%, the evolution is slower than a logarithmic. 

Before closing this Section, we would like to present 
a relationship between the critical value of the angular 
momentum L zc (that is the maximum value of the angu- 
lar momentum, for which, stars moving near the galactic 
plane are scattered, displaying chaotic motion, for a given 
value of the flatness parameter 6), when all other param- 
eters are kept constant. The results, which were found 
numerically, are given in Figures lla-b. For the elliptical 
galaxy of Fig. 11a the values of all the parameters are as 
in Fig. 4, while for the disk galaxy of Fig. lib the values 
of all the parameters are as in Fig. 6. Orbits starting in 
the upper part of the [b — L zc ] planes are regular, while 
orbits starting in the lower part of this plane including the 
line are chaotic. We see in Figure 12, that in both cases 
the relationship between the flatness parameter b and L zc 
in linear. An explanation of all the above relationships, 
will be given in Section 5. 

4. The S(g) dynamical spectrum 

During the last twenty years, dynamical spectra have 
been frequently used in Galactic Dynamics, as fast indi- 
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Figure 11: (a-b): Plot of the L zc vs. b, (a-left): for the elliptical galaxy and (b-right) for the disk galaxy. See text for details. 
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Figure 13: (a-b): (a-left): The S(g) spectrum for a regular orbit and (b-right): The corresponding L.C.E. 
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Figure 14: (a-b): The S(g) spectrum for an orbit producing two islands and (b-right): The two islands in the (r — p r ) phase plane of Fig. 6. 
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Figure 15: (a-b): The S(g) spectrum for an orbit producing three islands and (b-right): The three islands in the (r — p r ) phase plane of Fig. 
6. 



cators of regular or chaotic motion (see Contopoulos et Caranicolas & Vozikis 1999, Karanis & Caranicolas 2001, 
al. 1995, Contopoulos and Voglis 1996, Patsis et al. 1997, Caranicolas & Papadopoulos 2007, Caranicolas & Zotos 
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Figure 16: (a-b): Similar to Fig. 14a-b, but for an orbit producing nine small islands in the (r — p r ) phase plane of Fig. 4. 



2010). In this Section, we shall introduce and use, a new 
type of dynamical spectrum, which is the S(g) spectrum. 
We define the dynamical parameter g as 



9i = 



Pz 



(10) 



where ri,p r i,p z i are the successive values on the Poincare 
(r,p r ), z — 0,p z > 0, surface of section. We shall define 
the new dynamical spectrum of the parameter g, through 
its distribution function 



S(g) 



AiV(g) 
NAg ' 



(11) 



where AN(g) are the number of the parameters g in the 
interval (g 7 g + Ag), after N iterations. The reader can 
find more details and theoretical explanations, regarding 
the shapes and behaviors of the dynamical spectra, in the 
pioneer works of Froeschle et al. 1993, Voglis and Con- 
topoulos 1994 and Contopoulos ct al. 1997. By definition, 
the g parameter is based on a complicate combination of 
coordinates and momenta. The numerator of the parame- 
ter g ia a function of r and p r while, the denominator is a 
much more complex function of the same variables, since 
p z is related to r and p r through the energy integral (7). 
Due to the fact that the values of p Z i are the successive 
values on the Poincare surface of section, they can never 
obtain a zero value (always p z > 0), which can infinite the 
g parameter. There are many reasons for introducing the 
new S(g) spectrum. Three important reasons are that: (i) 
it can identify islandic motion of resonant orbits, (ii) it 
can help us to understand the evolution of sticky motion 
and (iii) it is much faster than the L.C.E, because it needs 
about 10 4 time units, while L.C.E needs about 10 5 time 
units, in order to obtain reliable results. 

Figure 13a shows the S(g) spectrum for an orbit, pro- 
ducing one of the invariant curves surrounding the central 
fixed invariant point shown in the phase plane of Fig. 6. 
Initial conditions are: ro = 4,p r o = 0,zo = 0, while the 
value of p z0 is always found from the energy integral (7). 
As expected, one can observe a well defined U type spec- 
trum, indicating regular motion. The spectrum shown in 
Fig. 13a is entirely in the g > domain. This means that 



the corresponding orbit is a typical orbit of the 1:1 reso- 
nance, forming an invariant curve surrounding the central 
invariant point of Fig. 6. The corresponding L.C.E for a 
time period of 10 5 time units, is given in Fig. 13b. Figure 
14a shows the S(g) spectrum producing two islands shown 
in Fig. 6, with initial conditions: r — l,p r o — 6.5, z = 0, 
while in Fig. 14b we see the two islands in the (r — p r ) 
phase plane. Figure 15a shows the S(g) spectrum, for an 
orbit producing three islands, shown in Fig. 6. Initial 
conditions are: ro = 9.34, p r o = 0, zq = 0. Here we ob- 
serve three U type spectra, which is as much the number 
of islands of Fig. 15b. Moreover, we see that the left 
and right spectra are exactly symmetrical about the 5 = 
axis, while the central spectrum lies on both sides of this 
axis. This indicates, that two of the islands are symmet- 
rical about the r - axis, while the third intersects the r - 
axis. In other words, we have the case of a quasi-periodic 
orbit, with a starting point on the r - axis. Figures 16a- 
b are similar to Figs. 15a-b, but for an orbit producing 
nine small islands shown in Fig. 4. Initial conditions are: 
r = 9.613, p ra = 0, z = 0. 

Thus, we came to the conclusion that the S(g) spec- 
trum is very useful, in order to identify islandic motion 
of resonant orbits. In this point, we should emphasize 
that, with some of previously used spectra (see Karanis & 
Caranicolas 2002 and references therein) we were able to 
detect islandic motion, but the number of spectra in some 
occasions was smaller than the number of islands, because 
symmetric islands produces identical spectra. Using the 
S(g) spectrum and particularly by adding the term —rip ri 
to its definition, we have managed to over through this 
drawback. 

Figure 17a shows the S(g) spectrum for a chaotic orbit. 
In this case, we see a complicated, asymmetric spectrum, 
with a lot of large and small peaks. This orbit belongs 
to the chaotic sea shown in Fig. 6. Initial conditions are: 
ro = 3.4, p r o — 22, zo = 0. The corresponding L.C.E for a 
time period of 10 5 time units, is given in Fig. 17b. 

In Fig. 16a, we have seen the S(g) spectrum for an 
orbit producing the set of nine small small islands, shown 
in Fig. 4. Figure 18a shows the S(g) spectrum for an orbit 
starting near the above described regular orbit. Initial 
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Figure 17: (a-b): (a-left): The S(g) spectrum for a chaotic orbit and (b-right): The corresponding L.C.E. 




Figure 18: (a-d): Evolution of the S(g) spectrum of a sticky orbit. Details are given in text. 



conditions are: r = 9.605, p r o — 0, z — 0. Here, we can 
see again nine spectra each one corresponding to an island. 
The basic difference between Fig. 16a and 18a is that in 
Fig. 18a we observe some additional small peaks. Those 
additional peaks, indicate sticky motion. It is well known, 
that in a dynamical system of two degrees of freedom, 
sticky orbits are those orbits, which stay for long time 
periods near the last K.A.M torus before they escape to 
surrounding chaotic sea (see Karanis & Caranicolas 2002). 

Let us now start to follow the time evolution of a sticky 
orbit, using the new S(g) spectrum. The sticky period is 
about 1750 time units. Figure 18b shows the S(g) spec- 
trum of the sticky orbit about 400 time units after the 
test particle has left the sticky region. Here the nine spec- 
tra have joined together producing three sets of spectra. 
Figure 18c shows the S(g) spectrum of the sticky orbit 
about 700 time units after the test particle has left the 
sticky region. Here the nine spectra have joined all to- 
gether producing a single spectrum. This spectrum has 



the characteristics of a chaotic spectrum and strongly in- 
dicates, that the test particle has left completely from the 
sticky region and has gone to the chaotic sea. The shape 
of the S(g) spectrum after the test particle has traveled 
for 20000 time units in the chaotic sea, is shown in Figure 
18d. 

5. A semi - theoretical approach 

In the present Section, we shall present some theoreti- 
cal arguments, together with elementary numerical calcu- 
lations, in order to explain the numerically found relation- 
ships, given in Figs. 9, 10 and 12. Since the potential of 
the nucleus (3) is integrable with spherical symmetry, the 
contribution of the total potential, to the chaotic regions 
observed in the (r—p r ) phase planes of the dynamical sys- 
tem, should derive mainly from the F zg force, which is the 
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Figure 19: A semi - theoretical plot of |-F Z9 | vs. b, for the elliptical 
and disk galaxy models. The values of all other parameters are given 
in the text. 



vertical force in the z - direction. The F zg force is 



bz (a + Vh 2 + bz 2 ) 



Vh 2 + bz 2 r 2 + (a + Vh 2 + bz 2 ) + c 2 



.(12) 



Figure 19 shows a plot of the \F Z 



b, for the ellip- 



tical and the disk galaxy models. As the scattering occurs 
near the nucleus, it must be at r < 1 and z < 1. The 
particular values of r and z are irrelevant. In this case, we 
choose: r = r = 0.2 and z = z = 0.1. We observe, that 
\F zg \, increases almost linearly with b, in both cases. One 
can see that, the pattern between \F zg \ and the flatness 
parameter b (shown in Fig. 19), is similar to those be- 
tween the maximal L.C.Es or the chaotic percentage A% 
and the flatness parameter b (shown in Figs. 9 and 10 
respectively), which have been obtained numerically. In 
both cases numerical outcomes suggest that, as the value 
of the flatness parameter b increases, the chaoticity of the 
dynamical system increases as well. Numerically found re- 
lationships given in Figs. 9 and 10 can be explained, if 
F zg force is the main contributor, to the chaotic regions 
observed in the (r — p r ) phase planes of the dynamical 
system. 

In order to explain the results shown in Fig. 12, we use 
an analysis similar to that used in Caranicolas & Innancn 
(1991). When a star approaches the dense nucleus, its 
momentum in the z - direction, changes according to the 
equation 



mAv z = (F zt ) At, 



(13) 



where m is the mass of the star, At is the duration of the 
encounter and < F zt > is the total average F z force. It 
was observed, that stars deflection into higher z proceeds 
in each case cumulatively, a little more with each succes- 
sive pass by the nucleus and not with a single dramatic 
encounter. It is assumed, that the star is scattered off the 
galactic plane after n (n > 1) encounters, when the total 
change in the momentum in the z - direction, is of order 



of the tangential velocity — L z /r. Therefore we have 

n n 

^Av zi ^(F zt )^^U- (14) 



m 



i=0 



i=0 



Let 



m 



i, 

Lz 
r 



F z t — F zn + F. 



zn ~ ± zg j 



(15) 



where F zn stands for the constant vertical force, coming 
from the dense nucleus, computed at a given point and for 
fixed values of the involved parameters. Then Eq. (14) 
becomes 



r 



F zg )T e . 



(16) 



Initially, if v z > the orbit is chaotic. Otherwise, 
v z <C v ^ and after successive encounters with the galactic 
center, v z will increase, until it reaches v$. Then the orbit 
will be chaotic. The time needed to reach the condition 
Vz > is T > T e , where T e is the sum of the encounter 
durations (15). Let Tf be the age of the galaxy, which 
corresponds here to the final time of numerical integra- 
tions. If T e > Tf, then some stars do not have enough 
time to be scattered. We conclude from Eq. (16), that if 
L z > L zc , where L zc is the angular momentum computed 
with T e = Tf, then there is still exist non-chaotic motions. 
Whereas for L z < L zc , all stars may have time to deploy 
chaotic motions. The star must go close to the nucleus, in 
order to be scattered. In this case, we may set r = r < 1 
and z — z < 1 in (16) and then we conclude in 



L zc ss k + Xb, 



where 



k = 
A = 



F zn T f r , 
F zg T f r 



(17) 



(18) 



are constants, at a given point P(ro,Zo). Relation (17) 
is very simplistic, but nevertheless reproduces the correct 
dependency on the flattening parameter b. Specifically, 
this expression explains the linear relationship between 
the critical value of the angular momentum L zc and the 
flattening parameter b, shown in Fig. 12. The reader can 
find similar analysis in Caranicaolas & Innanen (1991) and 
Karanis & Caranicolas (2001). The different slope in Fig. 
12, for the elliptical galaxy and disk galaxy model can be 
explained, because the values of k and A, of equation (17) 
are different for the two types of galaxies, depending on 
the remaining dynamical parameters, such as v ,a and h. 
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6. Discussion and conclusions 

In the present paper, we have studied the properties 
of motion, in a new galactic dynamical model. This new 
model, is a combination of the logarithmic and the Miyamoto 
- Nagai model (1975). The choice seems to be very suc- 
cessful for several reasons. The first reason is that, by 
changing the appropriate values of the parameters , one 
can obtain a variety of galactic models. Another reason, 
is that the dynamical system is simpler and displays a 
large variety of families of orbits, compared to other mod- 
els, such as the model proposed by Carlsberg and Innanen 
(1987), (see also Caranicolas & Innanen 1991, Caranicolas 
1997). 

It is well known, that galaxies evolve. Today, two main 
methods are generally accepted, in order to explain the 
results obtained from observations. In the first method, 
galaxies evolve with little influence from their environ- 
ment, while in the second, small galaxies assemble to form 
more complicated stellar systems (see Guiderdoni 2002 
and references therein). Observations show that, ellipti- 
cal galaxies are redder than disk galaxies, because they 
are made of old stars and there is no more gas in order to 
form young stars. On the contrary, disk galaxies are more 
blue, they contain younger stars and they have still gas to 
form them. Our galactic model, can be used in order to 
follow the galactic evolution from fiat disk systems to ellip- 
tical ones. This can be obtained, using a time-dependent 
form of potential (1), in which the parameters (vo,a, h) 
varies with time. Indeed, it seems that the chaotic region 
is more significant in disk galaxies. As a result, stars can 
escape from the galactic disk, so as the mass distribution 
changes and looks like more and more that of an elliptical 
galaxy. 

When studying Galactic Dynamics, it is invariably found 
that the stellar rotation velocity remains constant, or flat, 
with increasing distance away from the galactic center (see 
Catinella et al. 2006). This result, is highly counterintu- 
itive since, based on Newton's Law of gravity, the rota- 
tional velocity would steadily decrease for stars further 
away from the galactic center. By this particular argu- 
ment, the flat rotational curves, seem to imply that each 
galaxy must be surrounded by significant amounts of dark 
matter. It has been postulated and generally accepted, 
that the dark matter would have to be located in a mas- 
sive, halo enshrouding each galaxy (see Caranicolas & Zo- 
tos 2009, Caranicolas & Zotos 2010). The first real sur- 
prise in the study of dark matter lay in the outermost 
parts of galaxies, known as galaxy halos. Here, there is 
negligible luminosity, yet there are occasional orbiting gas 
clouds, which allow one to measure rotation velocities and 
distances. The rotation velocity is found, not to decrease 
with increasing distance from the galactic center, implying 
that the mass distribution of the galaxy, cannot be consid- 
ered, like the light distribution. The mass must continue to 
increase: since the rotation velocity satisfies O 2 = GM/r, 
where M is the mass within radius r, we infer that M in- 



creases proportionally to r. This rise appears to stop at 
about 50kpc, where halos appear to be truncated. We infer 
that the mass - luminosity ratio of the galaxy, including 
its disk halo, is about 5 times larger than estimated for 
the luminous inner region, or equal to about 50. 

The flat profile of the rotation curve shown in Fig. la, 
corresponding to our elliptical galaxy model, seems to be 
compatible with those of NGC 4261 (see Van der Marel, 
et al. 1990) and NGC 2778, NGC 3818, NGC 584, NGC 
3557, NGC 1553, NGC 5866, NGC 4594 and NGC 3115 
(see Davies et al. 1983). At the same time, the fiat pro- 
file of the rotation curve shown in Fig. lb, corresponding 
to our disk galaxy model, seems to be compatible with 
those of NGC 801, NGC 7331, NGC 224 (M31), NGC 
2590, NGC 2841, UMa: NGC 3953, UMa: NGC 3992, 
NGC 5907, NGC 1097 and NGC 4321 (see Brownstein 
and Moffat 2006). Here we must clarify, that when we 
say that the flat profiles of the theoretical rotation curves 
seem to be compatible with those from real galaxies, we 
mean that, if we superpose Fig. la-b with the rotation 
curves from the references, we will see that the theoreti- 
cal rotation curves of our present dynamical model pass 
through the observation data point, within the error bars. 
Additionally, both rotation curves shown in Figs, la-b, re- 
main flat with increasing distance from the galactic center. 
In both cases, we see that the asymptotic behavior of the 
theoretical rotation curves is very similar to the referenced 
rotation curves obtained by observation data. In Fig. la- 
fa we observe a bump for small values of r. This bump is 
caused, in both cases, by the dense spherical nucleus. Here 
we must note that, the above mentioned rotation curves of 
galaxies, were used in order to choose the appropriate val- 
ues of the parameters (vo,a, h) in both the elliptical and 
disk galaxy model. The chosen values of these parameters, 
secure positive value of density everywhere. Therefore, our 
new dynamical model, is realistic and can be used in order 
to describe real elliptical and disk galaxies. 

We have presented a description of the phase plane, for 
both the elliptical and the disk galaxy for low angular mo- 
mentum stars, that is L z < 10. The reason for doing that, 
is because in this case, we observe interesting families of 
regular orbits while, at the same time, a considerable part 
of the surface of section is covered by chaotic orbits. On 
the other hand, regular orbits have as conserved quantity 
some third integral, or quasi integral besides the energy 
and the z - component of the angular momentum. The 
reader can find more detailed information about the be- 
havior of low angular momentum stars in a series of papers 
(Caranicolas and Innanen 1991, Caranicolas 1997, Kara- 
nis & Caranicolas 2001, Caranicolas & Papadopoulos 2003, 
Caranicolas & Zotos 2010). Here, we must emphasize that 
the phase planes shown in Figs. 4 and 6, show that a test 
particle (that is a star with a negligible mass), with energy 
E = 486 in the elliptical galaxy model, shows the orbital 
behavior, displayed by the phase plane of Fig. 4, while a 
test particle with energy E = 1425 in a disk galaxy, shows 
the orbital behavior displayed by the phase plane of Fig. 
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6. We have chosen these values of energy, in each case, in 
order to cover a phase plane of about the same distance 
from the center of the galaxy. 

In order to keep things simple, we have kept some of the 
parameters of the system constant and studied the behav- 
ior of orbits varying only four basic parameters (vq, a, b, h). 
The parameter vq of our model, scales the rotation curve 
due to the potential V g , while the triplet (a,b,h) deter- 
mines the density distribution and the structure of the 
dynamical system (elliptical or disk galaxy). It was found 
that, both in the elliptical and disk galaxy model, there 
is a linear relationship between the L.C.E and the flatness 
parameter b. On the contrary, the relation between the 
chaotic percentage A% of the phase plane and the flatness 
parameter b is not linear but exponential. In both mod- 
els, the numerically found outcomes were explained, using 
some semi - theoretical arguments. In the same sense, we 
have explained the numerically obtained relationship be- 
tween the critical value of the angular momentum L zc and 
the flatness parameter b of the dynamical system. 

A very effective and reliable tool for the study of mo- 
tion, is the new S(g) spectrum, which allows us to detect 
islandic motion of resonant orbits, as it produces as much 
spectra as the total number of islands in the (r — p r ) sur- 
face of section. One more advantage of this spectrum, is 
that it can be used in order to calculate the sticky period 
of an orbit and also to follow its time evolution towards the 
chaotic sea. In this point, we must point out that all previ- 
ously used spectra cannot identify tiny islands of resonant 
orbits (see Fig. 16b). We have strong numerical evidence, 
that the new definition of the S(g) spectrum has the ability 
to detect islandic motion of complicated resonant orbits, 
which produce multiple islands in the Poincare surface of 
section. It is in our future plans, to use it in several and 
more complicated dynamical systems, with additional the- 
oretical work, in order to check if it gives always the correct 
number of spectra regarding the number of islands in the 
Poincare surface of section. 

Finally, we would like to point out, that the chaotic 
phenomena, are more significant, in the case where a mas- 
sive nucleus was present in the core of the galaxy and also 
the chaotic motion is observed only for the low angular mo- 
mentum stars (see Caranicolas and Papadopoulos 2003). 
Moreover, our numerical investigation reveals that, the de- 
gree of chaos in disk galaxies is much larger than compared 
to that in elliptical galaxies. This result is of particular 
interest, due to the fact that is generally accepted today, 
that massive objects inhabit the core of most galaxies (see 
Sprurzem 2003). 



very useful suggestions and comments, which improved the 
quality of the present paper. 
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